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Abstract. The construction of optimal template banks for matched-filtering 
searches is an example of the sphere covering problem. For parameter spaces 
with constant-coefficient metrics a (near-) optimal template bank is achieved by 
the j4* lattice, which is the best lattice-covering in dimensions n < 5, and is 
close to the best covering known for dimensions n < 16. Generally this provides 
a substantially more efficient covering than the simpler hyper-cubic lattice. We 
present an algorithm for generating lattice template banks for constant-coefficient 
metrics and we illustrate its implementation by generating A* template banks in 
n = 2, 3, 4 dimensions. 



1. Introduction 

The detection of gravitational waves (GWs) in tlie noisy data of detectors ideally 
requires the knowledge of the signal waveform, in order to coherently correlate the 
data with the expected signal by matched filtering. Depending on the type of 
astrophysical sources considered, however, one typically only knows a parametrized 
family of possible waveforms (or approximations thereof). The unknown parameters 
of these waveforms could be, for example, the frequency and sky-position of spinning 
neutron stars, or the masses and spins of inspiralling compact binary systems. 
Parameter spaces of such wide-parameter searches typically have between one and four 
dimensions, depending on computational constraints and the amount of astrophysical 
information available to constrain the search space a-priori. In the case of GWs from 
general binary systems, however, the number of dimensions of the parameter space 
could be as large as 17. 

Obviously one can only search a finite subset of points in this parameter space, and 
this subset constitutes the "template bank" . The templates must cover the parameter 
space, i.e., they must be placed densely enough that no signal in this space can lose 
more than a certain fraction of its power (called mismatch) at the closest template. 
However, coherently correlating the data with every template is computationally 
expensive and increases the expected number of statistical false-alarm candidates. An 
optimal template bank therefore consists of the smallest possible number of templates 
that still guarantees that the worst-case mismatch does not exceed a given limit. 

It was realized early on that a geometric approach is very useful to construct 
template banks, in particular the introduction of a parameter-space metric [3l [11] 
based on the mismatch. This provides a natural measure of distance in parameter 
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space and allows one to "correctly" place templates, in the sense that the maximal 
mismatch is not exceeded. Less attention, however, was devoted to the problem of 
optimally placing templates once the metric is known. Early works have sometimes 
used a hypcr-cubic template grid for illustrative purposes or the problem was 
incorrectly referred to as a "sphere packing problem" [121 E] • We see in the following 
that constructing an optimal template bank is an instance of the sphere covering 
problem, which is somewhat "dual" to the sphere packing problem. The full solution 
to the sphere covering problem in Euclidean space is only known in n = 2 dimensions, 
partial solutions (restricted to lattices) are known in n < 5 dimensions, while an 
optimal solution for higher dimensions is unknown (cf. [8l[T4]). The main motivation 
of the present work is to develop a general method for constructing efficient template 
banks in dimensions n < 17 by using the known results about Euclidean sphere 
covering. 

Previous related work on template banks includes studies to optimally cover non- 
flat two-dimensional parameter spaces arising in searches for GWs from inspiralling 
compact binary systems [H [2]. An interesting algorithm to construct a hexagonal 
{A2) template bank for 2D inspiral searches was described recently in [7]. Various 
codes exist within the LIGO Scientific Collaboration to generate hyper-cubic 
lattices (LALCreateFlatMeshO [10]), two-dimensional grids for non-constant metrics 
(LALCreateTwoDMeshO [lU]), and a tree-dimensional template bank based on the bcc- 
lattice (LALInspiralSpinBankO [10]), which is being used in a search for spinning 
binary inspirals on LIGO Data [l]. 



2. Template-based searches and parameter-space metric 

A wide class of searches for GWs can be characterized as template based, in the sense 
that one searches for signals belonging to a family of waveforms s{t; A), which depend 
on a vector of parameters {A}' = A\ The strain x{t) measured by a detector contains 
(usually dominating) noise n{t) in addition to possible weak GW signals s{t; As), i.e., 
x(t) = n(t) +s{t]Xs). One typically constructs a detection statistic, J-"(A; x) say, 
namely a scalar characterizing the probability of a signal with parameters A being 
present in the data x{t). Due to the random noise fluctuations n{t) in the data, the 
detection statistic is a random variable, and generally (assuming is unbiased) its 
expectation value ^{X; As) = E[J^{X; x)] has a (local) maximum at the location of the 
signal A = As, i.e., 

dT{X; X, 



dX 



-0. (1) 



Taylor-expanding the expected detection-statistic T in small offsets AA 
around the signal location As therefore reads as 

1 92J(A;As) 



.F(A;As) =.F(As;As 



2 dX' dXi 



AA* AA^ +0(AA^), (2) 

A, 



where the matrix of second derivatives of J- is negative definite. Here and in the 
following wc use automatic summation over repeated parameter indices i,j,.... We 
can introduce a mismatch m, which characterizes the fractional loss in the expected 
value of the detection statistic, at a parameter-space point A, with respect to the 
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signal location Ag, namely 

^ ^ ^(A.;A,) ^ ' 

Using the local expansion ([2]), we find 

m(A; A,) = .g„-(A,) AVAA^' + 0{AX^) , (4) 

where we defined the positive-definite metric tensor gij = —^didjJ-, and di = d/dX^. 
When searching a parameter space P(A% gij), we need to compute the detection 
statistic J-{x; A^) for a discrete set of templates Aj E P. Generally one can distinguish 
two different approaches to this problem: one is a random sampling of P using Markov- 
chain Monte-Carlo (MCMC) algorithms (e.g. see [6l [9]), and the other consists of 
constructing a template bank T = {A^} C P that covers the whole of P, in the sense 
that no point A g P exceeds a given maximal mismatch mmax to its closest template 
A^ e T, i.e., 

max min m(\: Ac) < mmax • (5) 

Here we focus on the construction of optimal template banks, namely those satisfying 
([5]) with the smallest possible number of templates A^. In the local metric 
approximation ([4]), each template A^ covers a region of parameter space, namely 

B^ = {XeF:g,j (A^) AA* AA^ < m,„ax , AA = A - A^ } , (6) 

which is a sphere of radius R — ^rrimax in the metric space P(A*, gij). We can therefore 
reformulate the definition of an optimal template bank as the set of (overlapping) 
spheres of covering radius R which cover the whole of P in the sense of (O with the 
smallest number of spheres. This is known as the sphere covering problem |8], not to 
be confused with the somewhat dual sphere packing problem, which seeks to pack the 
largest number of non-overlapping "hard" spheres into a given volume. 



3. The Euclidean sphere covering problem 

In this section we summarize the current status of the sphere covering problem as far 
as relevant for the construction of optimal template banks. There has been impressive 
progress in the study of the covering problem in recent years, e.g. see [5] for a general 
overview and [13] for a more recent update. Unfortunately, all of these studies are 
restricted to Euclidean spaces E" , while the metric parameter spaces of GW searches 
are often curved. In the following we will therefore make the assumption that P(A', gij) 
can be treated as at least approximately flat, or can be broken into smaller pieces that 
can be treated as nearly flat. If the curvature of the metric is too strong, i.e., if 
the curvature radius is comparable to the covering radius, it will be difficult to make 
use of the Euclidean covering problem, and a different approach such as a stochastic 
template bank or an MCMC sampling might be more effective. We further assume 
that we have found a coordinate system of P such that the metric components are 
(approximately) constant, i.e., 5ij(A) ~ constij, and for simplicity of notation we 
assume in this section (without loss of generality) that we have chosen coordinates 
in which the constant-coefficient metric is Cartesian, i.e., P = E"(a;*, Sij). 

A covering can consist of any arrangement of covering spheres, but currently all 
best coverings known are lattices, and we therefore restrict the following discussion to 
lattice coverings. 
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3.1. Basics on lattices 



An ?i-dimensional lattice A can be defined as a discrete set of points (forming an 
additive group) generated by 

l^^^eli^), with fez, (7) 

with summation over i = 1, . . . ,n, and where is a basis of the lattice. Note 

that it is sometimes convenient to express the n basis vectors in a higher-dimensional 
Euclidean space, i.e., generally we can have G E™ with m > n. When writing E" 
in the following we refer to the subspace of E™ containing the n-dimensional lattice 
A. The m X n matrix M°'i = l"^^-^ is called a generator matrix of the lattice, with the 
columns of M holding the m components of the n lattice basis vectors, so we can also 
write the lattice A as 

(8) 



The n X n matrix A = M"^ M is called the Gram matrix (where T denotes the 
transpose), which is symmetric and positive definite, and Aij 



(0 ■ '(i) 



/) !, /'^ 1° 



i.e., its coefficients are the mutual scalar products of lattice basis vectors. Each choice 




Figure 1. Hexagonal lattice (Aj) illustrating a 2-dimcnsional lattice covering. 
The shaded areas are different choices of fundamental regions for the lattice. FPl 
and FP2 are fundamental polytopcs ^ associated with different choices of lattice 
basis, WS is the Wigner-Seitz cell mil l, and R is the covering radius. 



of lattice basis defines a corresponding fundamental parallelotope (FP), namely 

FP ({Z(,)}) = {a; G E" : a: = 0' , < 0' < l} , (9) 

which is illustrated in figure [1] The FP is an example of a, fundamental region for the 
lattice, i.e., a building block containing exactly one lattice point, which fills the whole 
space E" when repeated. There are many different choices of basis and fundamental 
regions for the same lattice A, but they all have the same volume vol(A), given by 

vol(A) = VdetA, (10) 

and in the case where M is a square matrix wc also have vol(A) = detAf. One 
special choice of fundamental region is the nearest-neighbor region, often referred to 
as Dirichlet-Voronoi cell by mathematicians, and more commonly known as Wigner- 
Seitz cell or Brillouin zone by physicists, which is defined as 

WS(A) = {xe E" : ||a; - i/qH < \\x - z/flj, for aU 1^4 G A} , (11) 
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where j|a;j| ~ y/x ■ x is the standard Euchdean norm in E" . The vertices of the 
Wigner-Seitz ceh are by construction local maxima of the distance function of points 
in E" from the nearest grid point. The maximum distance of any point in E" to the 
nearest point of the lattice is called the covering radius R, which corresponds to the 
circumradius of WS, as seen in figure [T] 

Two lattices Ai and A2 with generator matrices Mi and M2 are equivalent if 
they can be transformed into one another by a rotation, reflection and change of 
scale, namely if the generator matrices satisfy 

Ah^cBMiU , (12) 

where c G R is a scale-factor, U is integer-valued detU = ±1, which accounts for 
different choices of basis vectors, and i? is a real orthogonal matrix, i.e., B = I. 
The associated Gram matrices are therefore related by 

A2=c^U^AiU, (13) 

and the fundamental volumes PH)) of the two lattices are 

vol(A2) =c"vol(Ai). (14) 

Let us consider as an example the 2-dimensional hexagonal lattice, illustrated in 
figure [TJ An obvious generator matrix is 

Mi=(l J/'V (15) 



V3/2 

corresponding to FPl in figure [1] However, sometimes it is more convenient to work 
with a generator matrix of the form 

/ 1 0\ 

M2=\ -I 1 , (16) 

V -1/ 

which has simpler coefficients, but uses a 3-dimensional representation of the 2- 
dimensional lattice with all lattice points lying in the plane x + y + z = 0. One 
can verify that these two representations are equivalent in the sense of (fT2)) . namely 
with 



c = V2, U=(l JJ ) , B=[ -I/V2 -1/V6 I . (17) 

V / 

Such a higher-dimensional representation of the generator matrix will be useful later 
for the description of the n-dimcnsional A* lattice. 

3.2. Known results on optimal sphere covering 

The efficiency of a sphere covering can be characterized by its thickness Q (sometimes 
also referred to as the covering density), which measures the fractional amount of 
overlap between the covering spheres, or equivalently the average number of spheres 
covering any point in E". This can be expressed as the ratio of the volume of one 
covering sphere to the volume of the fundamental region of the lattice, i.e.. 



I/V2 -1/V6 \ 
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where R is the covering radius and Vn is the volume of the unit-sphere in n dimensions, 
namely K = n''^^/T{n/2 + 1). We also use the normalized thickness or center density 
9, defined as 

0- ^ 



(19) 



which corresponds to the number of centers (i.e., templates) per unit volume in the 
case of i? = 1. Note that under a lattice transformation the covering radius R 
obviously scales as R2 ~ cRi, and we therefore see from that the thickness 
and (jl9|) is an invariant property of a lattice, i.e., 62 ~ Oi. The covering problem 
consists of finding the covering with the lowest center density 0. 

Kershner showed in 1939 (see [8]) that in n = 2 dimensions the most economical 
arrangement of circles covering the plane is the hexagonal lattice, which is equivalent 
to an A2 lattice. In dimensions n = 3,4,5 only the best lattice covering is known, 
and is given by A* in all three cases. In three dimensions, is also known as the 
body- centered- cubic (bee) lattice. Note that the best packing in n = 2 is also achieved 
by the hexagonal lattice, but for n = 3 the face-centered cubic (fee) lattice provides a 
denser packing than bee. In higher dimensions the best lattice coverings are currently 
still unknown, but the best coverings known can be found in see table 2 of |14| . and 
[15] provides for an up-to-date online version. As will become clearer in the following, 
the A* lattice, while no longer the "record holder" for most dimensions 5 < n < 17, 
is still close to the best currently known covering in all cases. In the following we will 
therefore mostly focus on the A* covering. The A* lattice has a center density of 

n{n + 2) ^ 



0(a:) = xA7tt|- 



12(n+l)J ' 

while for the hyper-cubic grid Z" the Wigner-Seitz cell is a unit hypercube, so 
vol(Z") = 1, and the covering radius R — y/n/2 is half the length of the diagonal. 
Therefore the center density is found as 6{17^) = 2^" n"/^, which is dramatically 
worse than A* in higher dimensions, as can be seen from the thickness ratio 

, , 0(Z«) 3"/2 /n+l\"^'n-oo 3"/2 

= 77777T = / -r — ~ "r^ ■ (21) 



9{Al) V«^V^ + 2/ ^/ne 



Table 1. Thickness ratio K(n) = and ■y(n) = 6»(bcst)/e(A* ) in 

dimensions n < 17. 



n 


2 


3 


4 


5 


6 


7 


8 


9 


K{n) 


1.3 


1.9 


2.8 


4.3 


6.8 


10.9 


17.7 


28.9 


7(n) 


1.0 


1.0 


1.0 


1.0 


0.97 


0.95 


0.86 


0.97 


n 


10 


11 


12 


13 


14 


15 


16 


17 


K{n) 


47.4 


78.2 


130 


216 


359 


601 


1007 


1692 


7(n) 


0.98 


0.88 


0.99 


0.86 


0.82 


0.86 


1.0 


0.68 



There is a theoretical lower limit on the thickness of any covering, the Coxeter-Few- 
Rogers (CFR) bound r„ (see [8]), i.e., 0n > Tn/Vn, where asymptotically r„ ~ njie^fe) 
for n — > 00. Figure [2] shows the normalized thickness as a function of dimension n for 
the A*^ and hyper-cubic Z" lattices, as well as the CFR bound and the best covering 
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Figure 2. Normalized covering tiiickness d as function of dimension n, for the 
hyper-cubic lattice (Z"), the A* lattice, the theoretical lower bound (CFR), and 
the best lattice coverings known. 



known. In table [T] we see that in dimensions n > 5 where A* has been superseded 
as the best covering [T5], the relative improvement 7(71) = 6'(best)/6'(A* ) in thickness 
is typically quite small. In particular, for 71 < 16 the improvement 7(71.) is typically 
less than 18%, while the advantage K{n) of A* compared to the hyper-cubic grid Z" 
grows large very rapidly, as seen in table [T] and figure [21 For practical simplicity we 
therefore propose to use A* as the covering lattice of choice. 



4. Lattice covering of template spaces 



Template counting 

The template spaces P(A*, g^) with constant-coefficient metrics gij only differ from 
the Cartesian case of the previous section by a simple coordinate-transformation. An 
infinitesimal parameter-space region d^X has a volume dV measured by the metric, 
namely dV = y^d"X, where g = detgij. The volume y of a finite region of parameter 
space is therefore 



V= dV = ^ / d"A, (22) 
Jp Jp 

where we used the fact that gij is a constant-coefficient metric. The number of 
templates dNp in dV is given by the inverse lattice volume, i.e.. 

Using the relation R = ^yriniax together with ((T8)l . p9)) . we find 

dN, = e m-^/^ dV =^ N, = e V9 / ^"A , (24) 

which generalizes template counting [Til IHl [1] to arbitrary lattices. 
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4-2. Practical implementation of lattice covering 



In this section we present a practical algorithm for generating lattices covering of 
given maximal mismatch mmax- The approach described here works for any lattice 
generator M, but in practice (cf. section 13. 2[) we will be most interested in the A* 
lattice. The generator for A* can be expressed (cf. [8]) as an {n + 1) x n matrix, 

/I 1 ... 1 ^ \ 

n-\-l 

-1 



\ 



1+ 
1 

n+1 
1 

n+1 



1 

n+1 



(25) 



n+1 / 



where the columns of M hold the n lattice basis vectors expressed in E"+^, i.e., 
i"^-) , with index conventions i, j = 1, . . . , n and a,b ^ 1, . . . , n + 1. The volume 



of the fundamental region and the covering radius for this generator are 

1 



vo1(A:) 



and R{A*J 



' n{n + 2) 
12(n+ 1) ' 



(26) 



which yields the (normalized) thickness d{A*^) given in ([^0]) . In order to generate such 
a lattice in a parameter space P(A*, g^ ), we need to express the generator M'^j in the 
A* coordinates, resulting in j , say, such that the lattice of templates is generated 

by ^ 

\i = J f , with I e Z" . (27) 

This coordinate transformation can be achieved in several steps: 

(i) Reduce the {n+1) x n matrix M°-j to a full rank generator, A'P j say, by expressing 
the lattice basis vectors in a Euclidean basis spanning the n-dimensional subspace 
E" of the lattice: a simple Gram-Schmidt procedure with respect to the Cartesian 
metric 6ab is used on the {'"jj} to generate an orthonormal basis satisfying 

h 



Sab e^j) e^^) — 5ij 



(28) 



The full-rank generator M*j is obtained from the components of the lattice 
vectors in this orthonormal basis, namely 

b 



1 * 



la 



(29) 



(ii) Translate the full-rank generator AP j from Cartesian coordinates into the 



coordinate system A* with metric 



For this we use another Gram-Schmidt 



orthonormalization with respect to the metric gij, with the lattice vectors {1^^^} 
as input to find an orthonormal basis {d^j^} satisfying 

9'j dli) c^(fe) = Sik ■ (30) 

This representation of an orthonormal basis in coordinates A* allows us to express 
the lattice vectors in these coordinates as 
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(iii) Scale the generator to the desired covering radius R = -^mmax, and with (|26p we 
find 



'l2(n + 1) 



If 



(32) 



/ n{n + 2) (j)' 

which is a generator (|27p for an A* template lattice with maximal mismatch 



This algorithm has been implemented in XLALFindCoveringGenerator () in LAL |10j . 

and some tests of this code are presented in the next section. 



4-3. Tests of the implementation 

In order to illustrate and test the implementation of this algorithm, we generate 
an A* lattice in dimensions n = 2,3,4, respectively, with a maximal mismatch of 
'TT-max = 0.04, i.e., a covering radius of i? = 0.2. For generality we use a non-Cartesian 
metric gij ^ Sij, as illustrated in the left panel of figure [3l We picked 100,000 points 




-0.6 -0.3 0.3 0.6 

X 



Figure 3. Left panel: Hexagonal (Aj) lattice covering in coordinates {x,y} with 
metric gij = [1, 0.4; 0.4, 0.5]. Right panel: Histogram of measured distances y/m, 
in a Monte-Carlo sampling of 100,000 points from an A* covering in n = 2, 3, 4 
dimensions, using non-Cartesian metrics gij. The nominal covering radius in all 
three cases was R - 



: 0.2. 



A G P(A*, gij) at random and computed their mismatch m (using the metric) to the 
nearest template , which is a way of measuring the maximal mismatch of a template 
bank. The distribution of measured mismatch-distances y/m is plotted in the right- 
hand panel of figure[31 and we see that the mismatches are bounded by ^Wmax = 0.2, 
satisfying We can also measure the (normalized) thickness 9 of the template bank, 
namely from the number of templates Np in the covered parameter space AA", we 
find using 



These measured values of the thickness arc found to agree to within 0.2% with the 
theoretical values (|20p in all three cases n = 2, 3, 4. The generated template banks in 
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this example have Np ~ 0(10"') templates, and the error can most likely be attributed 
to boundary effects. 

5. Discussion 

Possible applications of this algorithm for GW searches can be found in template- 
based searches, such as for inspiralling compact binary systems and for "continuous 
waves" , which in ground-based detectors refers mostly to signals from spinning neutron 
stars, and in the case of LISA includes white dwarf binaries, supermassive black hole 
binaries and extreme-mass ratio inspirals. The benefit of using this approach depends 
sensitively on the number of parameter-space dimensions, but can be estimated from 
table [1] at least in comparison to hypercubic grids. 

However, the applicability of the lattice covering algorithm presented here is 
restricted to explicitly fiat parameter spaces, which limits its usefulness to cases where 
we can find a coordinate system in which the parameter-space metric is (at least) 
approximately constant. The orbital metric approximation |13| for continuous GWs 
can be shown to be flat (work in progress), and would therefore be a natural case 
where this lattice covering could be used to greatest effect. One difficulty in this case, 
however, stems from that fact that the corresponding metric is found to be highly ill- 
conditioned, which results in the lattice-construction algorithm to fail due to numerical 
problems. One therefore needs to analytically "factor out" this near-degeneracy of the 
metric before this lattice-covering procedure can be safely applied. More work is also 
required to deal with non-trivial parameter-space boundaries, which complicates the 
n-dimensional filling algorithm. 
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